با استفاده از داده های زلزله ها در ایران و جهان به سوالات زیر پاسخ دهید.

۱. با استفاده از داده های historical_web_data_26112015.rds و استفاده از نمودار پراکنش سه بعدی بسته plotly نمودار طول، عرض و عمق زلزله ها را رسم نمایید. علاوه بر آن بزرگی هر نقطه را برابر بزرگی زمین لرزه قرار دهید. 

data = read_rds(paste(path, 'historical_web_data_26112015.rds', sep = '/'))
plot_ly(data, x= ~Latitude, y= ~Longitude, z= ~Depth, size= ~10*Magnitude) %>%
    add_markers

۲. پویانمایی سونامی های تاریخی را بر حسب شدت بر روی نقشه زمین رسم نمایید.(از داده زلزله های بزرگ استفاده نمایید.) 

disaster = read_delim(paste(path, 'disaster.txt', sep = '/'), "\t", escape_double = FALSE, trim_ws = TRUE)
disaster %>% select(I_D, Power = EQ_PRIMARY, COUNTRY, LOCATION_NAME, LATITUDE, LONGITUDE, REGION_CODE, DEATHS) -> 
    largeRep
mapWorld <- borders() # create a layer of borders
mp <- ggplot() +   mapWorld
mp + geom_point(data = largeRep, aes(x=LONGITUDE, y=LATITUDE, size=Power))

۳. نمودار چگالی دو بعدی زلزله های تاریخی ایران را رسم کنید.( از داده iran_earthquake.rds و لایه stat_density_2d استفاده نمایید). 

iran = read_rds(paste(path, 'iran_earthquake.rds', sep = '/'))
p = ggplot(iran, aes(x= Long, y= Lat))
p + stat_density_2d()

۴. احتمال اینکه در ایران در پنج سال آینده زلزله به بزرگی هفت ریشتر رخ دهد را محاسبه کنید. (از احتمال شرطی استفاده کنید.)

In this part data for earthquake disasters is used for isolating earthquakes with magnitude higher than 7 happened in Iran, then the time distance between these earthquakes have been calculated and used to find the distribution of earthquake time distance. The distribution is shown below: 

iran = disaster %>% select(COUNTRY, Power = EQ_PRIMARY, YEAR) %>% filter(COUNTRY == 'IRAN') %>% filter(YEAR >= 0) %>% 
    na.omit() %>% filter(Power >= 7)
iran$YEAR = as.numeric(iran$YEAR)
iran %>% mutate(dist = 0) %>% arrange(YEAR) -> iran
for(i in 0:nrow(iran)){
    if(i > 1){
        iran$dist[i] = iran$YEAR[i] - iran$YEAR[i-1]
    }
}
iran = iran[2:nrow(iran),]
q =  ggplot(data = iran) + stat_density(aes(x = dist), geom="line",position="identity")
q = print(q)

We see from the data that last earthquake was in 2017, so for meassuring the probability of an earthquake larger than 7 in the following five years can be calculated by integrating the above figure from time distance of 1 to 6 years, which is done below: 

q = q$data[[1]] 
q = q %>% select(x, y) %>% filter(x >= 1) %>% filter(x <= 6)
print(tail(cumtrapz(q$x, q$y), 1) * 100) 
##           [,1]
## [14,] 4.113475

So we see that the probability is 4.11%.

۵. بر اساس داده های زلزله های بزرگ ابتدا تعداد و متوسط کشته زلزله ها را بر حسب کشور استخراج نمایید. سپس نمودار گرمایی تعداد کشته ها را بر روی کره زمین رسم نمایید.(مانند مثال زیر!) 

largeRep %>% group_by(COUNTRY) %>% 
    summarize(deathMean = mean(DEATHS, na.rm = T), long = LONGITUDE[1], lat = LATITUDE[1], reg_code = REGION_CODE[1]) -> deathRep

mapData = get_data_from_map(download_map_data("custom/world"))
data = mapData %>% select(code = 'hc-a2', COUNTRY = name) 
data$COUNTRY = str_to_upper(data$COUNTRY)
data = full_join(data, deathRep %>% select(COUNTRY, deathMean))
data = data[1:213,]
data = data %>% select(code, deathMean) 
hcmap("custom/world", data = data, value = 'deathMean')

۶. با استفاده از داده لرزه های بزرگ و به وسیله طول، عرض، شدت، عمق مدلی برای پیش بینی تعداد کشته های زلزله بیابید.

Since the relationship between death and used features are highly non-linear, using linear model is not reasonable. So we are going to use another kind of regressor that works by ensembling a lot of trees and can model non-linear relationships. 

disaster %>% select(mag = EQ_PRIMARY, long = LONGITUDE, lat = LATITUDE, depth = FOCAL_DEPTH, death = DEATHS) %>% 
    na.omit()-> modRep
model = xgboost(data = as.matrix(modRep %>% select(1:4)), label = as.matrix(modRep %>% select(death)), 
                nrounds = 1500)
modRep$pred = predict(model, as.matrix(modRep %>% select(1:4)))
modRep$err = (modRep$death - modRep$pred)^2
rmse = sqrt(mean(modRep$err))
print(paste0("The root mean square of model is: ", rmse ))
## [1] "The root mean square of model is: 0.0413451777970055"

۷. با استفاده از داده worldwide.csv به چند سوال زیر پاسخ دهید. تحقیق کنید آیا می توان از پیش لرزه، زلزله اصلی را پیش بینی کرد؟

The foreshock happens in a range between 1 day to several years before the occurence of the main quake. So in many cases it is impossible to use foreshock for prediction of main quake. Some great earthquakes such as 1950 India-China earthquake with a magnitude of 8 didn’t have any foreshocks at all.

۸. گزاره " آیا شدت زلزله به عمق آن بستگی دارد" را تحقیق کنید؟ (طبیعتا از آزمون فرض باید استفاده کنید.) 

world = read_csv(paste(path, 'worldwide.csv', sep = '/'))
world %>% select(depth, mag) %>% na.omit() -> magDepthRep
cor.test(magDepthRep$depth, magDepthRep$mag)
## 
##  Pearson's product-moment correlation
## 
## data:  magDepthRep$depth and magDepthRep$mag
## t = 50.03, df = 60629, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.1914584 0.2067469
## sample estimates:
##       cor 
## 0.1991147

۹. میانگین سالانه زلزله ها را بر حسب کشور به دست آورید. آیا میتوان دلیلی در تایید یا رد تئوری هارپ ارائه کرد. 

world$year = substr(world$time, 1, 4)
world %>% select(year, place) -> rep
rep$country = str_match(rep$place, "(?<=, ).*$")
rep %>% group_by(country, year) %>% summarize(count = n()) -> rep

The countries that claim to be the victom of HAARP project are Iran and Pakistan. The plot for number of earthquakes per year in these 2 countries are shown below: 

world$year = substr(world$time, 1, 4)
world %>% select(year, place) -> rep
rep$country = str_match(rep$place, "(?<=, ).*$")
rep %>% group_by(country, year) %>% summarize(count = n()) -> rep
hchart(rep %>% filter(country == 'Iran'), hcaes(x = year, y = count), type = 'column') %>% 
    hc_yAxis(title = list(text = 'Number of Earthquakes')) %>% hc_title(text= 'Iran') %>% 
    hc_xAxis(title = list(text = 'Year'))
hchart(rep %>% filter(country == 'Pakistan'), hcaes(x = year, y = count), type = 'column') %>% 
    hc_yAxis(title = list(text = 'Number of Earthquakes')) %>% hc_title(text= 'Pakistan') %>% 
    hc_xAxis(title = list(text = 'Year'))

If the HAARP project was really affecting these countries earthquake numbers, we should see an increasing pattern in the figures, which we don’t. So it’s just a conspiracy theory.

۱۰. سه حقیقت جالب در مورد زلزله بیابید.

I. In Japan mythology a giant catfish called Namzu is responsible for earthquakes. In Hindu mythology, Earth is held in place by eight gigantic elephants, all balanced on the back of turtle, which itself stands on the coils of a snake, if any of these animals move an earthquake will happen.

  1. The most powerful earthquake ever recorded was in Chile in 1960 and had a magnitude of 9.5.

  1. The deadliest earthquake known hit Shansi, China, 1556 and killed an estimated 830000 people.